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We use molecular dynamics simulation results on viscous binary Lennard-Jones mixtures to ex- 
amine the correlation between the potential energy and the virial. In accord with a recent proposal 
[U. R. Pedersen et al. Phys. Rev. Lett. 100, 015701 (2008)], the fluctuations in the two quantities 
are found to be strongly correlated, exhibiting a proportionality constant, F, numerically equal to 
one-third the slope of an inverse power law approximation to the intermolecular potential function. 
The correlation is stronger at higher densities, where interatomic separations are in the range where 
the inverse power law approximation is more accurate. These same liquids conform to thermody- 
namic scaling of their dynamics, with the scaling exponent equal to P. Thus, the properties of strong 
correlation between energy and pressure and thermodynamic scaling both reflect the ability of an 
inverse power law representation of the potential to capture interesting features of the dynamics of 
dense, highly viscous liquids. 

PACS numbers: 61.43.Fs, 61.20.Lc, 64.70.Pf, 61.20.Ja 



I. INTRODUCTION 

The myriad phenomena accompanying the supercool- 
ing of hquids continue to intrigue condensed matter sci- 
entists. On approaching the glass transition, the equi- 
librium, metastable liquid is comprised of fast and slow 
moving molecules, which interchange their roles at times 
on the order of the structural relaxation. This makes the 
many-body dynamics inherently cooperative and thus 
too complex to be solved by a first-principles treatment. 
Various models have been developed to describe and in- 
terpret the dynamic properties, but none rise to the level 
of a predictive theory. Nevertheless, progress can be 
made by combining experiments that capture some defin- 
ing characteristic of the dynamics with computer simu- 
lations able to provide a microscopic interpretation. 

One effective simplification in describing the dynamics 
of glass-forming liquids is provided by the experimen- 
tal observation that structural relaxation times Tq, and 
transport coefficients, encompassing a range of thermo- 
dynamic conditions, superpose when plotted as a func- 
tion of p^/T, where T is temperature, p density, and 7 
a material constant [H, [H, @] . This thermodynamic scal- 
ing of Tq, has been demonstrated for dozens of molecular 
liquids and polymers [4] and extends with good accuracy 
from the high temperature Arrhenius region, through the 
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dynamic crossover, down to the glass transition [s^. The 
justification for the scaling law 

r„ ~ HpVT) (1) 

draws from early work of Hoover and coworkers 0], 
who investigated the properties of fiuids having a pair- 
wise additive intermolecular potential described by an 
inverse power law (IPL) 

u{r) = e{a/rr (2) 

in which m is a constant, r the intermolecular separation, 
and e and a have respective dimensions of energy and 
length. Generally, the intermolecular potential for liquids 
is represented as a pair- wise additive interaction, with the 
steep repulsive part written as an exponential function 
or an IPL. An advantage of the IPL is that all reduced 
thermodynamic and dynamic properties of such a system 
can be expressed in terms of the variable p^^^/T 
Simulations of vitrifying liquids have often employed an 
IPL potential [§, [lO, [ill, [l2[- Recent theoretical models 
of energy landscapes, relevant for description of the glass 
transition, are also based on Eq. ([2]) but with the addition 
of a density-dependent constant to account for the long- 
range attractive forces [mil [Hill. 

At high densities the liq uid structure is determined 
mainly by the repulsions [17] , suggesting that if the repul- 
sive potential were accurately described by an IPL, the 
local dynamics would depend only on p^/^/T, with the 
empirical scaling exponent 7 identified with m/3. Sup- 
port for this idea comes from recent simulations using an 
m-6 Lennard-Jones (LJ) intermolecular potential 
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Relaxation times and diffusion constants from these sim- 
ulations conform to Eq. ([1]) b ut with a scaling exponent 
7 that is not equal to m/sfisl. fl9l.[2Q|. For polymers this 
is due in some measure to the effect of the intramolecular 
part of the potential f^lj. More generally, a systematic 
study of LJ particles [20] showed that 7 follows closely 
the steepness of the repulsive core. In the latter sys- 
tem, however, the scaling exponent 7 is always larger 
than m/3, due to contributions of the attractive term in 
Eq. da]) [H, m, [23. 7 is equal to the slope of an IPL fit- 
ted over the typical distances of closest approach between 
particles probed in the highly viscous regime. These in- 
teratomic separations are defined by values of r between 
the first nonzero value of the radial distribution function 
and the position of the half- height of the first peak [2Q| |. 

In addition to giving rise to the scaling of the dy- 
namics [Eq. ([1])], adequacy of the IPL approximation for 
highly viscous liquids suggests the existence of a corre- 
lation between equilibrium fluctuations of the configu- 
rational parts of energy and pressure; i.e., between the 
potential energy U and the virial W. This expecta- 
tion follows from the exact correlation between fluctu- 
ations in U and W for an IPL [2^. Recently, Dyre and 
coworkers showed from simulations that the potential en- 
ergy and the virial strongly correlate for various non- 
associated materials, with correlation coefficients exceed- 
ing 0.9 [13, Hsf. Liquids exhibiting both thermodynamic 
scaling and correlation between U and W may exhibit 
other interesting properties, such as sufficiency of a sin- 
gle parameter to describe their temperature-dependent 
viscoelasticity [26|, [l^ . In this work we assess the corre- 
lation of U and W in viscous LJ liquids having different 
repulsive interactions; i.e., different m in Eq. (|3]). These 
are the same LJ liquids previously shown to conform to 
thermodynamic scaling of their diffusion coefficients [lO] • 
Thus, from our results we can test the conjecture of Ped- 
ersen et al. [25^ that liquids whose dynamics follow the 
thermodynamic scaling are strongly correlating and vice 
versa. 



II. MODELS 

The models considered in this work consist of additive, 
equimolar binary mixtures of 500 particles interacting 
through the m-6 LJ potential 
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where a,/? = 1,2 are indexes of species. The mix- 
tures studied have a size ratio crii/<j22 = 0.64, masses 
mi/ 7X12 = 1-0, and a unique energy scale 612/611 = 
622/611 = 1.0. We refer to them herein as AMLJ. As 
in previous work [20], we consider different repulsive in- 
teractions: m = 8, 12, 24, and 36. We also study two 
other LJ liquids having m = 12: one proposed by Wahn- 
strom [Isl, which we denote as WAHN, and the non- 
additive mixture of Kob and Andersen [2^]. This latter. 
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FIG. 1: Normalized instantaneous fluctuations of the po- 
tential energy {U — (U)) / y^{AlP) (circles) and virial {W — 
{W))/y^{AW^) (squares) for the AMLJ model with m = 36 
at the indicated pressure and temperature. The ordinate of 
the upper curves has been shifted 5 units for clarity. 



called herein BMLJ, is among the most widely used in 
molecular dynamics simulations of the glass transition. 
In the followi ng reduced LJ units will be used, assuming 
CTii, 611 and ^/mlall/ell as units of distance, energy and 
time respectively. All samples are cooled isobarically at 
pressures P = 5, 10, and 20 using Berendsen thermo- 
stat and barostat during equilibration. The production 
runs are then performed in the NVE ensemble using the 
Velocity- Verlet algorithm. Further details on the simula- 
tions of AMLJ models can be found in Refs. 

We remark that the actual pressure varies with m even 
though the numerical values of P are the same, because 
the depth of the potential increases with m. A more 
appropriate set of reduction parameters for the poten- 
tial in Eq. (|4]) is given by the position r* and the width 
ix* = 1^11 (r*) of the minima of unir). In terms of these 
parameters, the effective pressure is P* = 
where P is expressed in reduced LJ units. Considering 
the variation of u* and a* for different m, the pressure 
P for m = 36, for example, should be increased by about 
a factor of 2.5 to match that of m = 12. 



III. RESULTS 

Figure [T] is a plot of the instantaneous values of the 
normalized fluctuations in the virial and the potential 
energy for the AMLJ liquid with m = 36 for two se- 
lected conditions of T and P. As discussed below, this 
material exhibits the weakest correlation of U and W 
within the class of systems and range of state conditions 
investigated herein; nevertheless, a relationship between 
the two quantities is evident. The degree of correlation 
can be assessed by a variety of means. We calculate the 
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FIG. 2: Fluctuations of the virial versus those of the potential energy in the AMLJ model with m = 8, 12, 24, and 36. For 
each m, results for all pressures P — 5, 10, and 20 and all temperatures studied are shown. Colors in the online version 
indicate datasets corresponding to different state points. The corresponding range of densities is 1.04 < p < 2.06 for m = 8, 
0.87 < p< 1.78 for m = 12, 0.93 < p< 1.70 for m = 24, and 0.77 < p< 1.71 for m = 36. The solid lines represents the least 
squares linear fits with the correlation coefficient R indicated. 



Pearson correlation coefficient 

V((Af/)2){(AW-)2) 

from linear regression of AU = U — (U) and AW = 
W — (W)^ with the large number of data points, typi- 
cally 2 X 10^ per sample, enhancing statistical reliability. 
The determination of R is carried out simultaneously for 
all state points over which thermodynamic scaling was 
observed in Ref. [20]; that is, for each pressure (P = 5, 
10, and 20) at temperatures corresponding to normal liq- 
uid conditions down to the slow-dynamics regime. 

In Fig. [2] the fluctuations in the virial are plotted ver- 
sus those in the potential energy for AMLJ liquids with 
m = 8, 12, 24, and 36. The obtained R are listed in 
Table [H from which two observations can be made: The 
correlation coefficients are all close to unity, indicating 
substantial correlation. This concurs with the results of 



Bailey et al. [23| for simulations of one-component LJ 
liquids with m = 12 and of BMLJ particles. Second, 
the magnitude of R decreases from ca. 0.97 to 0.83 with 
increasing magnitude of the repulsive exponent. This 
ostensibly suggests that liquids with steeper intermolec- 
ular repulsions are less strongly correlating. However, as 
discussed in Section [III the actual pressure regime ex- 
plored by the particles in the simulation changes with m. 
Thus, larger m corresponds to smaller effective P and p. 
Since the correlation improves at higher density [23], re- 
flecting the greater accuracy of the IPL approximation, 
the smaller R for larger m can be ascribed to the dif- 
ferences in thermodynamic conditions. We demonstrate 
this in Fig. [3l showing results for m = 36 calculated for 
P = 30 and 50. At these higher pressures the correla- 
tion coefficient is similar to the value of R obtained for 
m = 12 at P = 10 and 20. This confirms that the poorer 
correlation (smaller R) at higher m is due to the lower 
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TABLE I: Summary of exponents from U-W correlations (F) 
and from thermodynamic scaling of diffusion (7) for various 
LJ liquids. The correlation coefficient R obtained from linear 
regression of W{U) data is also included. Statistical uncer- 
tainties on r correspond to one standard deviation in the 
"mean" case, and to the error associated to linear regression 
in the "global" case. 

m/3 7 r (global) R T (mean) 

AMLJ36 12 13.4 ± 0.2 12.21 ± 0.02 0.827 12.1 ± 0.9 
AMLJ24 8 9.1 ± 0.1 9.07 =b 0.01 0.872 9.0 dz 0.2 
AMLJ12 4 5.0 ± 0.1 5.095 ± 0.004 0.934 5.10 ± 0.17 
AMLJ8 2.7 3.5 ± 0.1 3.667 ± 0.002 0.965 3.67 ± 0.13 
BMLJ 4 5.0 ± 0.1 5.087 ± 0.003 0.943 5.10 ± 0.15 
WAHN 4 5.0 ± 0.1 5.052 ± 0.003 0.978 5.16 ± 0.19 



FIG. 3: Fluctuations of the virial versus those of the potential 
energy for the AMLJ model with m = 36 at pressures P = 30 
and 50 and ah T studied (yielding 1.20 < p < 1.75). These 
correspond to pressures P* (see Section |Il} comparable to the 
values for m = 12 at P = 10 and 20. 
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FIG. 4: Least squares slope F calculated along individual 
isobars shown as a function of density p = p{T,P). Results 
are shown for the AMLJ model with m = 8 (0.18 < T < 2.0), 
12 (0.41 < T < 3.0), 24 (0.75 < T < 4.0), and 36 (0.87 < T < 
6.0). The horizontal dashed lines denote the corresponding 
values of the scaling exponent 7 (from Ref. [20| V 



effective pressure, and thus larger mean nearest neighbor 
distances, rather than to the steepness of repulsion per 
se. 

For an exact IPL ^ = m/3 [23] and regression of 
W{U) yields the value of m/3 as the slope of the fitted 
line. More generally, while the fact that R ^ 1 affirms 
correlation between U and VK, this is not a proof of pro- 
portionality of the two quantities 0|. However, it can be 



observed that the scatter in the plots of Fig. [2] is normally 
distributed (random scatter without systematic trends), 
thereby justifying an interpretation of the slope, F, of 
the fitted lines as a measure of m/3. Results for all sim- 
ulations, including the BMLJ and WAHN, are listed in 
Table [H where it can be seen that there is good corre- 
spondence between the slope F and the scaling exponent 
7; that is, strongly correlating liquids conform to thermo- 
dynamic scaling, thus confirming the results of Ref. [33^ . 

We can examine the relationship between pressure and 
energy fluctuations in more detail by evaluating the cor- 
relation for each state point individually. Figure [H shows 
F(T, P), obtained from linear regression, as a function of 
p = p(T, P). Results for AMLJ systems with m = 8, 12, 
24, and 36 are shown in the figure, along with the scal- 
ing exponents (indicated by the dashed line) obtained 
from superpositioning of the diffusion constants for these 
liquids. Interestingly, at fixed m, F(T, P) essentially col- 
lapse onto a single horizontal curve when plotted versus 
p, with only some variation at lower T and P. These 
changes in F(T, P) reflect the fact that the IPL approxi- 
mation depends weakly on the state point, with the fluc- 
tuations of F(T, P) for different state points mostly dic- 
tated by density variations (again excepting m = 36 for 
the reasons discussed above). For each m, the mean of 
F(T, P) over all state conditions is equivalent within the 
error to the F obtained from a global fit of AW versus 
AU (see Table III). 

One final comment concerns the link between the F 
from the U-W correlations and the thermodynamic scal- 
ing exponent 7. In Ref. poj we found that a single value 
of 7 gave excellent superpositioning of the dynamic data. 
This means that any change in "local" 7 with T and P 
must be small (this issue was examined quantitatively 
for simulated liquid silica data in Ref. [s^])- Even for the 
more poorly correlating liquid (m = 36), the estimated 
uncertainty in 7 is only about 15% (±0.2). On the other 
hand, the slopes describing U-W correlations display a 
somewhat wider variation upon changing state param- 
eters, especially at low density and temperature. This 
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may be related to the sensitivity of pressure-energy fluc- 
tuations to the shape of the pair potential for distances r 
around and beyond the flrst peak in the radial distribu- 
tion function (e.g., departures from the IPL form) [23 |. 
Further investigations on the role of the attractive tail 
of the potential on U-W correlations may be required to 
clarify this point. 

IV. CONCLUSIONS 

In summary, we have shown that for viscous m-6 LJ 
liquids conforming to thermodynamic scaling of their dy- 
namics, there is a strong correlation between the virial 
and the potential energy. This property is maintained for 
systems ranging from relatively soft particles (m = 8) 
to those approaching the hard sphere hmit (m = 36). 
The correlation deteriorates for lower densities, as the 
range of the fluctuations extends to large r, for which 
the IPL approximation breaks down. The slopes ob- 



tained from linear regression of the virial vs. potential 
energy data are in good agreement with the scaling ex- 
ponents yielding superpositioning of dynamic data, sup- 
porting the conjecture that pressure-energy correlations 
and thermodynamic scaling have a common origin in the 
IPL approximation of the interaction potential. Our re- 
sults are in accord with the recent work of Dyre and 
coworkers [23, 24, 25, 3^ on LJ particles and suggest the 
utility of the IPL approximation in describing essential 
features of the dynamics of dense, highly viscous liquids. 

Acknowledgments 

Computational resources were obtained through a 
grant within the agreement between the University of Tri- 
este and CINECA (Italy). D. C. acknowledges financial 
support by the Austrian Science Fund (FWF) (Project 
number: P19890-N16). The work at NRL was supported 
by the Office of Naval Research. 



[1] R. Casalini and C. M. Roland, Phys. Rev. E 69, 062501 

(2004) . 

[2] C. Alba-Simionesco, A. Cailliaux, A. Alegria, and G. Tar- 
jus, Europhys. Lett. 68, 58 (2004). 

[3] C. Dreyfus, A. Le Grand, J. Gapinski, W. Steffen, and 
A. Patkowski, Eur. J. Phys. 42, 309 (2004). 

[4] C. M. Roland, S. Hensel-Bielowka, M. Paluch, and 
R. Casalini, Rep. Prog. Phys. 68, 1405 (2005). 

[5] R. Casalini and C. M. Roland, Phys. Rev. B 71, 014210 

(2005) . 

[6] W. G. Hoover and M. Ross, Contemp. Phys. 12, 339 
(1971). 

[7] W. G. Hoover, S. G. Gray, and K. W. Johnson, J. Chem. 

Phys. 55, 1128 (1971). 
[8] Y. Hiwatari, H. Matsuda, T. Ogawa, N. Ogita, and 

A. Ueda, Progr. Theor. Phys. 52, 1105 (1974). 
[9] B. Bernu, J. P. Hansen, Y. Hiwatari, and G. Pastore, 

Phys. Rev. A 36, 4891 (1987). 
[10] T. Grigera, A. Cavagna, L Giardina, and G. Parisi, Phys. 

Rev. Lett. 88, 055502 (2002). 
[11] C. De Michele, F. Sciortino, and A. Coniglio, J. Phys.: 

Condens. Matter 16, L489 (2004). 
[12] L. A. Fernandez, V. Martin-Mayor, and P. Verrocchio, 

Phys. Rev. Lett. 98, 085702 (2007). 
[13] P. G. Debenedetti, F. H. Stillinger, T. M. Truskett, and 

C. J. Roberts, J. Phys. Chem. B 103, 7390 (1999). 
[14] M. S. Shell, P. G. Debenedetti, E. La Nave, and 

F. Sciortino, J. Chem. Phys. 118, 8821 (2003). 
[15] R. J. Speedy, J. Phys.: Condens. Matter 15, S1243 

(2003). 

[16] M. S. Shell and P. G. Debenedetti, Phys. Rev. E 69, 
051102 (2004). 

[17] J. D. Weeks, D. Chandler, and H. C. Andersen, J. Chem. 

Phys. 54, 5237 (1971). 
[18] J. Budzien, J. D. McCoy, and D. B. Adolf, J. Chem. 



Phys. 121, 10291 (2004). 
[19] G. Tsolou, V. A. Harmandaris, and V. G. Mavrantzas, 

J. Chem. Phys. 124, 084906 (2006). 
[20] D. Coslovich and C. M. Roland, J. Phys. Chem. B 112, 

1329 (2008). 

[21] C. M. Roland, S. Bair, and R. Casalini, J. Chem. Phys. 

125, 124508 (2006). 
[22] D. Ben-Amotz and G. Stell, J. Chem. Phys. 119, 10777 

(2003). 

[23] N. P. Bailey, U. R. Pedersen, N. Gnan, T. B. Schr0der, 

and J. C. Dyre, J. Chem. Phys. 129, 184507 (2008). 
[24] N. P. Bailey, U. R. Pedersen, N. Gnan, T. B. Schr0der, 

and J. C. Dyre, J. Chem. Phys. 129, 184508 (2008). 
[25] U. R. Pedersen, N. P. Bailey, T. B. Schr0der, and J. C. 

Dyre, Phys. Rev. Lett. 100, 015701 (2008). 
[26] N. L. Ellegaard, T. Christensen, P. V. Christiansen, N. B. 

Olsen, U. R. Pedersen, T. B. Schr0der, and J. C. Dyre, 

J. Chem. Phys. 126, 074502 (2007). 
[27] N. P. Bailey, T. Christensen, B. Jakobsen, K. Niss, N. B. 

Olsen, U. R. Pedersen, T. B. Schr0der, and J. C. Dyre, 

J. Phys.: Condens. Matter 20, 244113 (2008). 
[28] G. Wahnstrom, Phys. Rev. A 44, 3752 (1991). 
[29] W. Kob and H. C. Andersen, Phys. Rev. E 51, 4626 

(1995). 

[30] D. Coslovich and G. Pastore, J. Chem. Phys. 127, 124504 
(2007). 

[31] D. Coslovich and G. Pastore, J. Chem. Phys. 127, 124505 
(2007). 

[32] F. J. Anscombe, Am. Stat. 27, 17 (1973). 

[33] T. B. Schr0der, U. R. Pedersen, and J. C. Dyre, 

arXiv:0803.2199 (2008). 
[34] A. Le Grand, C. Dreyfus, C. Bousquet, and R. M. Pick, 

Phys. Rev. E 75, 061203 (2007). 



